#ifndef MODEL_FACTORIZATION_MACHINE_CC
#define MODEL_FACTORIZATION_MACHINE_CC

#include "model_fm.hpp"
#include "def.hpp"
#include "utils.hpp"
#include "log.hpp"
#include <cmath>
#include <functional>
#include <numeric>
#include <algorithm>

using namespace std;

namespace Puza{
	PuzaModelFM::~PuzaModelFM(){ }

	PuzaModelFM::PuzaModelFM(){
		// initialize all latent factors
		PuzaLogger::Instance()->PutString("(Factorization Machines):Init latent factors.");
		this->theta.resize(PuzaDEF::Instance()->INPUT_ROW_TABLE.size());
		for(unsigned i=0; i < this->theta.size(); i++){
			this->theta[i].resize(PuzaDEF::Instance()->FACTOR);
			for(int k=0; k < PuzaDEF::Instance()->FACTOR; k ++){
				this->theta[i][k] = gaussian_random_variable(0, PuzaDEF::Instance()->INIT_STDEV);
			}
		}

		// initialize all weights
		PuzaLogger::Instance()->PutString("(Factorization Machines):Init all other weights.");
		int J = PuzaDEF::Instance()->FEATURE_MAPPING.size();
		this->weights.resize(J);
		for(int k=0; k < J; k++)
			this->weights[k] = gaussian_random_variable(0, PuzaDEF::Instance()->INIT_STDEV);

		// set global bias
		PuzaLogger::Instance()->PutString("(Factorization Machines):Init global bias.");
		this->weight_0 = 0.0;
	}

	double PuzaModelFM::Squared_Loss(double y_bar, double y){
		return (y_bar - y) * (y_bar -y);
	}

	double PuzaModelFM::Logistic_Loss(double y_bar, double y){
		double result = log(1 + exp(-y * y_bar));
		return result;
	}

	double PuzaModelFM::Huber_Loss(double y_bar, double y){
		double temp = y_bar * y;
		double result = 0.0;
		if(temp > 0){
			result =  0.5 * pow(max(0.0, 1 - temp),2);
		}
		else{
		  result = 0.5 - temp;
		}
		return result;
	}

	double PuzaModelFM::PointwiseLossFunction(double y_bar, double y){
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_SQUARE)
			return this->Squared_Loss(y_bar, y);
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_LOGISTIC)
			return this->Logistic_Loss(y_bar,y);
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_HUBER)
			return this->Huber_Loss(y_bar,y);
		return 0.0;
	}

	double PuzaModelFM::Hypothesis(PuzaRow r){
		double result = this->weight_0;
		// calculate linear regression part
		for(unsigned int i=0; i < r.features.size(); i ++){
			 int feature_id = r.features[i].feature_id;
			 double feature_value = r.features[i].feature_value;
			 result = result + this->weights[feature_id] * feature_value;
		}
		// calculate pairwise latent factors
		for(unsigned int i=0; i < r.features.size(); i++){
			for(unsigned int j= i + 1; j < r.features.size(); j++){
				int feature_id_1 = r.features[i].feature_id;
				int feature_id_2 = r.features[j].feature_id;
				double interaction_result = r.features[j].feature_value * r.features[j].feature_value;
				interaction_result = interaction_result * inner_product(this->theta[feature_id_1].begin(), this->theta[feature_id_1].end(), this->theta[feature_id_2].begin(), 0);
				result = result + interaction_result;
			}
		}
		return result;
	}

	void PuzaModelFM::Gradients(double temp, PuzaRow r){
		double lambda = PuzaDEF::Instance()->LEARNING_RATE;
		double reg = PuzaDEF::Instance()->REG;
		// gradient descent for global bias
		double weight_0_gradient =  temp + 2 * reg * this->weight_0;
		this->weight_0 = this->weight_0 - lambda * weight_0_gradient;
		// gradient descent for feature coefficients, including latent factors
		for(unsigned i=0; i < r.features.size(); i++){
			double feature_id = r.features[i].feature_id;
			double weight_k_gradient = temp * r.features[i].feature_value + 2 * reg * this->weights[feature_id];
			this->weights[feature_id] = this->weights[feature_id] - lambda * weight_k_gradient;
			// latent factors
			for(int k=0; k < PuzaDEF::Instance()->FACTOR; k++){
				double temp_result = 0.0;
				for(unsigned j=0; j < r.features.size(); j++){
				if(j != i){
					temp_result = temp_result + this->theta[j][k] * r.features[j].feature_value;
				}
			}
			double v_i_k_gradient = temp * (r.features[i].feature_value * temp_result) + 2 * reg * this->theta[i][k];
			this->theta[i][k] = this->theta[i][k] - lambda * v_i_k_gradient;
			}
		}
	}

	void PuzaModelFM::Squared_Gradients(double y_bar, PuzaRow r){
		double temp = 2 * (y_bar - r.response);
		this->Gradients(temp, r);
	}

	void PuzaModelFM::Logistic_Gradients(double y_bar, PuzaRow r){
		double temp = (sigmoid_function(r.response * y_bar)-1.0) * r.response;
		this->Gradients(temp, r);
	}

	void PuzaModelFM::Huber_Gradients(double y_bar, PuzaRow r){
		double temp = y_bar * r.response;
		if((temp < 1)&&(temp > 0)){
			temp = (temp - 1) * r.response;
		}
		if(temp >= 1){
			temp = 0.0;
		}
		if(temp <= 0){
			temp = - r.response;
		}
		this->Gradients(temp, r);
	}

	void PuzaModelFM::PointwiseSGD(double y_bar, PuzaRow r){
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_SQUARE)
			this->Squared_Gradients(y_bar, r);
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_LOGISTIC)
			this->Logistic_Gradients(y_bar, r);
		if(PuzaDEF::Instance()->CURRENT_POINT_LOSS == LOSS_HUBER)
			this->Huber_Gradients(y_bar, r);
	}
};
#endif
